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1.  INTRODUCTION 

Recently,  Jong-Sen  Lee  [1]  proposed  an  algorithm  for  contrast  enhancement 
and  white  noise  filtering  of  a  digital  image  based  on  local  statistics. 

Lee’s  algorithm  is  nonrecursive  and  requires  large  computational  effort  in 
the  calculation  of  local  means  and  variances  for  individual  pixels. 

In  what  follows, we  present  a  recursive  algorithm  for  the  calculation 
of  the  above  means  and  variances,  which  reduces  considerably  the  number  of 
computations  required  to  process  given  sets  of  pixels. 

2.  BAS  1C  ASSUMPT  IONS 

Let  x ( i , j )  denote  the  intensity  (grey  level)  of  the  pixel  (i,j).  Lee 
[1]  makes  the  following  assumptions: 

(a)  x(i,j)  is  a  realization  of  a  real-valued  random  variable  X ( T  ,  j ) 
whose  mean  and  variance  we  denote  by  ^(i,j)  and  <j2(  i , j )  ; 

(b)  the  observation  z(i,j)  of  x(i,j)  consists  of  x(i,j)  corrupted  by 
white  additive-  or/and  multiplicative  noise,  the  mean  and  variance  of  the 
noise  being  assumed  known; 

(c)  observations  z(k,l)  of  plxels(k^)  located  in  a  (2n  +  1)  x  (2m  +  1) 

rectangular  window  centered  around  a  given  pixel  (i,j)  are  independent  and 

* , 

all  come  from  the  same  scalar  distribution.  Hence  estimates  J  (i,j)  and 

*2  2  ^ 

aZ  (J»j)  of  Mz( i ,j)  andcz(i,j)  are  obtained  by  local  averaging  over  those 

pixels,  i.e. 

*  j  i+n  j+m 

M7  (l,J)  =  (2n+l)(2m+l)  ,  ?  E  z(k,t),  (]) 

L  k= i -n  T = j -m 

*2  ]  *+n  j+m 

°Z  °’J)  =  72n+i )  ( 2m+l  j"  ,  ?  ^  (z(k^)  -  J  (i,j))2.  (2) 

k='-n  t= j-m  Z 

For  each  pixel,  Lee  computes  uz,(i,j)  and  az2(i,j)  by  means  of  (1) 
and  (2)  and  then  uses  these  values  in  obtaining  a  minimum  mean  square 
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estimate  x(t,j)  of  x(i,j),  as  explained  (for  the  additive  noise  case) 
in  sect  ion  4. 

Actually,  the  observation  z(i,j),  should  not  be  included  in  the 
above  estimates,  since  it  is  used  later  on  (see  section  4)  in  the  esti¬ 
mation  of  x(i,j).  For  this  reason,  in  the  present  paper,  we  replace  (1) 
and  (2)  by  (la)  and  (2a)  below: 

*  j  i +n  j  +m 

MZ<'’J>  =  (2n+l)(2m+l)-l[  ?  2  z(k»-t)  -  z(i,j))  (la) 

k=  i -n  j -m 


'*2  1  *  +n  j  +m  „ 

0rz(I,J)  =  (2n+l )  (2m+l )- 1  <  *  S  (z(k.O  -  u  (i,j))2  -  (r(i,J)-y  (i,j))2)  (2 
4  k=  i  -n  ),-]  -m  Z  7 


There  are  some  other  issues  that  Lee's  approach  raises  which  we  will 

not  discuss  here.  The  only  objective  of  the  present  paper  is  to  develop 

a  recursive  algorithm  for  the  calculation  of  J  (i,j)  and  *  2(i,j)  defined 

Z  Z 

by  (la)  and  (2a). 


3.  RECURSIVE  ESTIMATION  OF  LOCAL  MEAN  AND  VARIANCE 

Let  the  domain  of  the  image  be  N  X  M.  Our  recursion  applies  to  a 
concentric  (N-2n)  x  (M-2m)  rectangular  subset  S  of  the  image  (see  Fig. 

1).  Thus  the  recursion  begins  at  the  pixel  (n+1 ,  m+1)  from  where  it 
moves  horizontally  rightward  up  to  the  pixel  (n+1,  M-m)#  Then  it  shifts 
to  the  next  horizontal  line  moving  rightward  from  the  pixel  (n+2,  m+1)  to  pixel 
(n+2,  M-m).  This  procedure  is  repeated  up  to  the  (N-n)th  horizontal  line. 


3.1.  Recursion  for  the  Local  Mean 


Let  us  define  the  following  entities 


C  O  ,j)  = 


(2n+l ) (2m+l ) - 1 


z(  i  ,j) 


(3) 


4 

i+n 

pO  J)  "  E  £d<.j)  .  '  >  n  +  1  ,  (4) 

k=  i  -  n 

7  ( i » j )  =  p(»,j)  -  £(i-n,j),  i  >  n  +  1  ,  (5) 


jj  (  ' .  j)  =  (a  (  i ,  j)  +  £(  i  ,j) 

z  z 


i  >  n  +  I,  j  >  m  +  I  ,  (6) 

Substituting  (3)  and  (4)  in  (la),  we  obtain 
.  j-hn 

U  (i,j)  =  E  p(i,/,)  -  £(i,j)  ■  (7) 

.  Z  ^=j -m 

A  detailed  consideration  of  the  above  variables  vis-a-vis  formula 
(7)  leads  to  our  proposing  the  following  recursion  for  the  generation 

A 

of  u  ( i  >j)  • 

Z 


(a)  Initiali zat i on 

(I)  For  j=! , . . . ,2m+! ,  calculate: 
Step  1:  CO.j),  i=] , . . .  ,2n+l , 

Step  2 :  p(n+l ,j) . 

( I  I)  Calculate: 

Step  3:  u  (n+1 ,  m+1) 

Z 
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(b)  Remaining  elements  of  (n+l)th  row 

For  j  =  m+2 , . . . ,M-m,  calculate 

Step  4:  q( i ,j+m) ,  i=l,...,2n+l  ; 

Step  5:  p  (n+1 ,j+m) 

Step  6:  pr  (n+1  ,j-l) 

Z 

Step  7:  |a^(n+l,j)  =  ^  (n+l,j-l)  +  p(n+l,j+m)  -  r(n+l,j) 

t  h 

B-  q  row,  q  ■=  n+2  , . . .  ,N-n 
(a)  Initial i zat i on 

I.  For  j=l , . . . , 2m+l ,  calculate: 
iStep  8:  Q (q+n ,  j ) 

Step  9:  p'(q-l.j) 

Step  10:  p  (q ,  j )  -jr(q-IJ)  +  e(q+n,j) 


II.  Calcul ate: 

Step  11:  ^  (q,  m+1) 

Z 

(b)  Remaining  elements  of  qth  row 


For 

j 

=  nrt-2  , . . 

calculate: 

Step 

12: 

q  ( q+n ,  j  +m) 

Step 

13: 

P 

(q-1 , j+m) 

Step 

14: 

p(q,j+m)  = 

IT  (q-1, j+m)  +  £ (q+n, j+m) 

Step 

15: 

pr 

7 

(q,j-D 

Step 

16: 

L. 

a 

U 

7 

(q.j)  = 

U  (q’j"1)  +  p(q , j+m)  - 

6 


3.2.  Recursion  for  the  Local  Variance 


Equation  (2a)  can  be  rewritten  in  the  form 
*  ^  j  i+n  j+m 

°7  <ilJ)  '  1  (in+lXlSl)-!  <  ,  ?  T-  *  (k.t)  -  22(i.j)) 

L  k=  i  -n  £=j-m 

.2 

-  |J  O.j)]  .  (8) 


This  shows  that  a  recursion  for  ^  (i,j)  can  be  formulated  which 

Z 

is  almost  identical  to  that  for  the  local  mean  with  z(k,/.)  replaced 
2 

by  z  (M).  In  fact,  let 


11  (l,j)  (2  n+1)  (2m+l)- 

2 

n  z  o,j) 

(9) 

i+n 

Y  O.J)  =  S  T|(k,j)  , 

i  >  n+1 

(10) 

k=  i  -n 

Y  ( '  J)  =  y  (  '  ,j)  -  ti  0- 

'  n , j ) ,  i  >  n+1 

(ID 

2  ^  2 

o  ( i » j )  =  a  ( i .  j )  +  rj(  i » j ) 
Z  Z 


y( i , j -m) 


r 


*  2 

+  (i  (  i  ,j)  ,  i  >  n+1  ,  j  >  m+1 
Z 

Note  that  (8)  may  be  rewritten  as 


,  2  i4™ 

<*  O.j)  =  £  y(i,T) 

z  -f  =j  -m 


-  J  2  (i,j) 

Z 


02) 


03) 
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We  may  now  proceed  to  state  the  algorithm  for  the  recursive  generation 

A 

of  a  (  i ,  j)  as  fol lows. 

Z 

A.  n  +  row 

(a)  Initial izat ion 

1.  For  j  =  1 , . . . , 2m+ 1 ,  calculate: 


Step  1  : 

T)(  1  *  j  )  »  i  -  li...  ,2n+l ; 

Step  2: 

y(n+l ,  j ) 

II. 

Calculate: 

Step  3: 

*2 

0  (n+1,  m+1)  (by  means  of  formula  (13),  where  in 

7 

> 

L. 

A 

this  formula  (n+1,  mfl)  is  obtained  from  Step 

Z 

in  the  algorithm  of  the  preceding  subsection). 

*  t"  h 

(b)  Remaining  elements  of  (n+1;  row 

For 

j  =  m+2,...,M-m,  calculate: 

Step  4: 

T]  ( i ,  j+m)  ,  i  =  l ....  ,  2n  +  1  ; 

Step  5: 

y  (n+1 ,  j+m) 

Step  6: 

J  2  (n+1,  j-1) 

7 

Step  7: 

L 

02  (n+1 ,  j  )  =  o’2  (n+1,  j-1)  +  y(n+l  ,  j+m) 

Z  Z 

-  T]  (n+1  ,j)  -  J  2  ( n+1  ,j) 
Z 


(where  £  (n+l,j)  is  obtained  from  Step  7  in 

Z 

the  algorithm  of  the  preceding  subsection). 
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(a)  Initial izat i on 

I.  For  j  =  1 , . . . , 2m  +  1  ,  calculate 

Step  8:  r|  (q  +  n,j) 

Step  9:  y  (q  -  1 ,j) 

Step  10:  Y  (q,j)  =  y(q-l,j)  +  ^  (q  +  n,j) 


I  I .  Calculate 
*  2 

Step  11:  (q,  mfl)  (by  means  Df  f:ormu]a  (]3)  where,  in 

this  formula,  Q  (q,  m+1)  is  obtained  from  Step  11 
.  in  the  algorithm  of  the  preceding  section). 


(b)  Remaining  elements  of  gth 


For  j  =  mf2,. . . ,M-m,  calculate: 


Step  12:  -q  (q+n,  j+m) 

Step  13:  y  (q-1 •  j+m) 

Step  14:  y  (q>j+m)  =  y  (q-U  j+m)  +  n  (q+n  ,  j+m)  > 

Step  15:  a2(q,j-l) 

Z 

Step  16:  a  (q,j)  =  a  (q, j-l )  +  y(q,j+m)  -  ■qCq.j)  -  J  2  (q,J) 

1  1  1 

.  *  2 

(where  (q,j)  is  obtained  from  Step  16  in  the 

algorithm  of  the  preceding  section). 
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4.  CONCLUSION  OF  THE  ALGORITHM  AND  FINAL  REMARKS 

The  local  mean  and  variance  obtained  as  above  can  now  be  used  in 
Lee  s  [I]  procedure.  To  show  this,  we  briefly  discuss  the  case  of 
filtering  in  the  presence  of  white  additive  noise.  For  the  cases  of 
multiplicative  noise  and  joint  additive  and  multiplicative  noise,  the 
reader  is  referred  to  [1], 

In  the  additive  noise  case  then,  we  write 

z( ' i j)  =  x(  i ,  j )  +  w( i , j)  ,  ( ]4) 


where  w(i,j)  is  a  realization  of  W(i,j),  which  is  independent  of 

X(i,j).  with  zero  mean  and  variance  ffy2  (Note  also  assumption  (c)  in 
sect  ion  2), 

Assuming  that  the  "a-priori"  mean  and  variance  of  X(i,j)  are 
respect ively 


and 


wx°'J) 


A 

n  ('.j) 
Z 


a  <  i,j)  =  C7 

X  z 


(  i  >  j  ) 


(15) 


(16) 


Lee  deduces  the  a-posteriori  estimate  x(?,j)  of  x(i,j)  to  be 
x('.j)  =  J^('.j)  +  K(i,j)(z(i,j)  -  Mx  (  »  ,  j  ) ) 


where 


K(i,j)  = 


aY  ( i ,  j ) 


2  2 
aX  +  aW 


(17) 


(18) 
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So,  In  this  case,  we  would  simply  substitute  y  (i,j)  and 
2  Z 

$  (l,j),  obtained  in  the  preceding  sub-sections,  into  (15)  and  (16), 

and  then  proceed  to  (17)  and  (18). 

At  this  point,  it  is  worthwhile  to  compare  the  computational 
effort  in  executing  the  recursive  (given  above)  and  non-recursive 
schemes  for  the  estimation  of  local  means  and  variances.  For  this 
purpose,  the  number  of  scaling  (multiplication  by  a  constant),  addition, 
and  multiplication  operations  for  the  (n+!)th  line  and  for  the  qth  line, 
where  q  >  n  +  1,  are  listed  in  Table  1  for  each  case.  Typically,  m  =  n  =  3 
and  M  =  N  =  256.  Table  2  shows  the  numbers  of  operations  listed  in 
Table  1  for  these  values  of  m,  n,  M  and  N. 

Typically,  in  the  calculation  of  all  the  local  means  and  variances 
associated  with  the  qth  line,  where  5^qs  253,  the  number  of  scaling 
operations  as  well  as  the  number  of  multiplications  for  the  nonrecursive 
and  recursive  algorithms  are  nearly  the  same.  However,  (for  the  qth  line) 
the  recursive  algorithm  requires  only  about  1/8  of  the  number  of  additions 
required  by  the  nonrecursive  one  in  the  calculation  of  the  local  means 
and  about  1/6  of  the  number  of  additions  needed  by  the  nonrecursive 
algorithm  in  the  evaluation  of  local  variances. 

From  the  above  we  conclude  that  in  processing  a  picture  based  on 
local  statistics,  one  could  combine  recursion  with  parallel  processing 
either  to  simplify  the  parallel  processor  without  sacrifice  in  speed,  or 
to  speed  up  the  overall  processing. 


REFERENCE 

kidi  Leoi’r,n0n9“Sen’  "Digital  lma9e  Processing  by  Use  of  Local  Statistics' 

NRL  Report  8192,  Washington,  D.C.,  March  1,  1978  ?o  appear  in  the 

Transactions  on  Pattern  Analysis  and  Machine  Intelligence. 


TABLE  1 


1  t 


NUMBERS  OF  ELEMENTARY  OPERATIONS 
FOR  PROCESSING  A  SINGLE  LINE 

(Estimation  applies  to  pixels  (i,m+l)  to  (i,M-m)  for  the  ith  line) 


Type  of  Operation 


Sea  I fng 


Addition  * 


Local  Means 


No.  of  Elementary  Operations 


L  i  ne 


(n+l)1"^  1  ine 

th  . 
q  1  me 

(n+2  <  q  <  M-m) 
(n+l)t^1  line 

th  . 

P  line 

(n+2  <  q  <  M-m) 


Non-Recurs i ve 
Algor i thm 
(J.S.  Lee) 

(We  use  (la)  and  (8) 
instead  of  (1)  and  (2)) 

M  -  2m 
M  -  2m 


((2n+l )  (2rrrt-l  )  -2)  • 
(M-2m) 

(  (2n+l ) (2m+l )-2) . 
<'M-2m) 


Recurs i ve 
Algor i thm 
(Present  Paper) 


(2n+l)M 

M 

2(n+2)M- 

3(2m+l) 

6M-3(2m+l ) 


Local  Var i ances 


Seal ing 


Add i t ion 


Mul t ipl i cat  ion 


(n+1) line 
th 

q  line 
(n+2  <  q  <  M-m) 

(n+1) ^  line 


th  . 
q  line 

(n+2  <  q  <  M-m) 


(n+1) ^  line 
qth  line 
(n+2  <  q  <  M-m) 


M-2m 

M-2m 


E ( 2 n+ 1 ) (2m+l ) - 1  ] . 
(M-2m) 


E ( 2  n+ 1 ) (2m+!) -1 j . 
(M-2m) 

(2n+l)M+M-2m 

M+M-2m 


( 2  n+ 1 ) M 
M 

(2n+6)M- 
5(2rrrf  ])  +  ] 


8M-5(2m+l)  +  1 


2(n+l )M-2m 
2(M-m) 
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TABLE  2 

VALUES  OF  NUMBERS  OF  ELEMENTARY  OPERATIONS 
Indicated  in  Table  1,  with  m=n=3,  M=N=256 

No.  of  Elementary  Operations 


Type  of  Operation 

Line 

Non-Recurs ive 
Algorithm 
(J.S.  Lee) 

Recurs ive 
Algor ithm 
(Present  Paper) 

Loca  1  s 

Means 

Seal  ing 

,  th  . 

4  1  me 

250 

1,792 

th  .  . 
q  line 

(5  <  q  <  253) 

250 

256 

Add  i  t ion 

,  th  . 

4  line 

11,750 

2,539 

* 

th  .  . 
q  1  ine 

(5  <  q  <  253) 

11,750 

1  ,515 

Local  Var iances 


Seal ing 

,  th  .  . 

4  1  ine 

250 

1,792 

,th  .  . 
q  1  ine 

(5  <  q  <  253) 

250 

256 

Add i t ion 

.  th 

4  line 

12,000 

3,038 

i 

th  . 
q  line 

(5  <  q  <  253) 

12,000 

2,014 

\ 

Mul t ipl i cat  ion 

,  th  .  . 

4  1  ine 

2,042 

2,042 

th  .  . 
q  1  me 

(5  <  q  <  253) 

506 

506 

Mul t ipl i cat  ion 


